/*-------------------<-- Start of Description-->---------------------\
| Binomial Equivalence Test; |
|---------------------<-- End of Description-->----------------------|
|--------------------------------------------------------------------|
|-----------<-- Start of Files or Arguements Needed-->---------------|
| Arguments: |
| This fuction take either variables from a dataset or direct |
| numbers; |
| indata: input data set; |
| ne: experiment group; |
| pe: success rate in experiment group; |
| se: count of success in experiment group; |
| nc: control group; |
| pc: success rate in control group; |
| sc: count of success in control group; |
| delta: = critical difference: |
| H0: (Control group - Experimental Group)=delta; |
| default is 0; |
| tail: = One Sided Test or Two Sided Test; |
| better:= the bigger the proportion the better, or the smaller |
| the better; BIG or SMALL: default is BIG; |
| Note: this option is required for One sided test; |
| test: HA (stands for Hauck-Anderson), BW (for Blackwelder), FM |
| (for Farrington Manning Procedure); |
| default is HA; |
| Note: the HA procedure One-Sided test will be the same as |
| BW procedure. |
| outdata: output dataset; |
|------------<-- End of Files or Arguements Needed-->----------------|
|--------------------------------------------------------------------|
|------------------<-- Start of Files Created-->---------------------|
| Example: %bnmleqv(ne=50, se=40, nc=60, sc=49,delta=0.10, test=FM, |
| tail=two, outdata=one); |
| Example 2: |
| data one; |
| success_e=70; |
| n_e=100; |
| success_c=49; |
| n_c=60; |
| run; |
| %bnmleqv(indata=one, ne=n_e, se=success_e, nc=n_c, sc=success_c, |
| test=fm, better=small, alpha=0.05, delta=0.1,tail=one); |
| Usage: %bnmleqv(indata=, ne=, pe=, se=, nc=, pc=, sc=, tail=one, |
| outdata=); |
\------------------<-- Start of Files Created-->--------------------*/
%macro bnmleqv(indata=, ne=, pe=, se=, nc=, pc=, sc=, delta=0, tail=one,
better=big,alpha=0.05,test=, outdata=);
/*--------------------------------------------\
| Author: Duo Zhou, Brian Johnson; |
| Created: 2-1-2002 7:02pm; |
| Purpose: Binomial Equivalence Test; |
\--------------------------------------------*/
%local _nnumeric_ _pnumeric_ _snumeric_ _consistence_ _i_ _grp_1 _grp_2 _linesize_ _tmplast_;
%let _tmplast_=&syslast;
%local _n1 _p1 _s1 _n2 _p2 _s2;
%let _grp_1=experimental;
%let _grp_2=control;
%let _linesize_ = %SYSFUNC(GETOPTION(linesize));
%let _n1=≠ %let _p1=&pe; %let _s1=&se; %let _n2=&nc; %let _p2=&pc; %let _s2=≻
%if (%quote(%upcase(&tail)) = ONE) or (%quote(%upcase(&tail)) = 1) %then %do;
%if (not %index(%quote(%upcase(&better)), BIG)) and (not %index(%quote(%upcase(&better)), SMALL)) %then %do;
%put ==> Alert! In order to setup the hypothesis, I need to know if the bigger the proportion;
%put +++ the better, or the smaller the proportions the better. Please specify a valid;
%put +++ option "BIGGER" or "SMALLER"!;
%end;
%end;
data _tmp1_;
retain _hyp0 hyp1 _diffSE _tvalue _pvalue;
%if &indata ne %then %do;
set &indata;
%end;
length _hyp0 _hyp1 $200.;
format _n1 _n2 _s1 _s2 Best10. default= 11.8;
%do _i_=1 %to 2;
%if (%quote(&&_n&_i_) ne) and (%quote(&&_p&_i_) ne) and (%quote(&&_s&_i_) ne) %then %do;
%if (%sysevalf(%sysevalf(&&_n&_i_*&&_p&_i_)-&&_s&_i_)>=0.5) %then %do;
%put ==> Alert! Inconsistency: (&&_n&_i_*&&_p&_i_ <> &&_s&_i_).;
%goto finish;
%end;
_n%trim(%left&_i_=&&_n&_i_;
_p&_i_=&&_p&_i_;
_s&_i_=&&_s&_i_;
%end;
%else %if (%quote(&&_n&_i_) ne) and (%quote(&&_p&_i_) ne) and (%quote(&&_s&_i_) eq) %then %do;
_n&_i_=&&_n&_i_;
_p&_i_=&&_p&_i_;
_s&_i_=(_n&_i_)*(_p&_i_);
%end;
%else %if (%quote(&&_n&_i_) ne) and (%quote(&&_p&_i_) eq) and (%quote(&&_s1) ne) %then %do;
_n&_i_=&&_n&_i_;
_s&_i_=&&_s&_i_;
_p&_i_=(_s&_i_)/(_n&_i_);
%end;
%else %if (%quote(&&_n&_i_) eq) and (%quote(&&_p&_i_) ne) and (%quote(&&_s1) ne) %then %do;
_p&_i_=&&_p&_i_;
_s&_i_=&&_s&_i_;
_n&_i_=(_s&_i_)/(_p&_i_);
%end;
%else %do;
%put ==> Alert! I will need at least 2 of the three parameters:_n&_i_,_p&_i_ and _s&_i_;
%put +++ in order to proceed.
%goto finish;
%end;
%end;
_p=(_s1+_s2)/(_n1+_n2);
%if (%quote(&test) eq) or (%quote(%upcase(&test)) eq %quote(BW)) or
%index(%quote(%upcase(&test)),%quote(BLACKWELDER)) or
(%index(%quote(%upcase(&test)),HA)) or
(%index(%quote(%upcase(&test)),%quote(HAUCK)) and
%index(%quote(%upcase(&test)),%quote(ANDERSON)))%then %do;
%if ((%index(%quote(%upcase(&test)),HA)) or (%index(%quote(%upcase(&test)),%quote(HAUCK))
and %index(%quote(%upcase(&test)),%quote(ANDERSON)))) and
((%index(%quote(%upcase(&tail)),TWO)) or (%quote(&tail) eq 2)) %then %do;
_diffSE=SQRT(_p1*(1-_p1)/(_n1-1)+_p2*(1-_p2)/(_n2-1));
%end;
%else %do;
%if (%quote(&delta) eq) or (%quote(&delta) eq %quote(0)) %then %do;
%let delta=0;
_diffSE=SQRT(_p*(1-_p)*(1/_n1+1/_n2));
%end;
%else %do;
_diffSE=SQRT(_p1*(1-_p1)/_n1+_p2*(1-_p2)/_n2);
%end;
%end;
_c=1/(2*min(_n1, _n2));
%if (%index(%quote(%upcase(&tail)),TWO)) or (%quote(&tail) eq 2) %then %do;
_hyp0="H0: |pc-pe|<=%trim(%left(&delta))"; _hyp1="Ha: |pc-pe| > %trim(%left(&delta))";
%if (%index(%quote(%upcase(&test)),HA)) or (%index(%quote(%upcase(&test)),%quote(HAUCK))
and %index(%quote(%upcase(&test)),%quote(ANDERSON))) %then %do;
_ucl = (_p1-_p2) + probit(1 - &alpha) * _diffSE + _c;
_lcl = (_p1-_p2) - probit(1 - &alpha) * _diffSE - _c;
if (1-probnorm((_p1-_p2+&delta-_c)/_diffSE)) >=(1-probnorm(abs((_p1-_p2-&delta+_c)/_diffSE)))
then _tvalue = (_p1-_p2+&delta-_c)/_diffSE;
else _tvalue=(_p1-_p2-&delta+_c)/_diffSE;
_pvalue=max((1-probnorm(abs(_p1-_p2+&delta-_c)/_diffSE)),(1-probnorm(abs((_p1-_p2-&delta+_c)/_diffSE))));
%end;
%else %do;
_lcl = (_p1-_p2) - probit(1 - &alpha/2) * _diffSE;
_ucl = (_p1-_p2) + probit(1 - &alpha/2) * _diffSE;
_tvalue = (_p1-_p2+&delta)/_diffSE;
_pvalue=(1-probnorm(abs((_p1-_p2+&delta)/_diffSE)))*2;
%end;
keep _lcl _ucl;
%end;
%else %if (%index(%quote(%upcase(&tail)),ONE)) or (%quote(&tail) eq 1) %then %do;
%if (%index(%quote(%upcase(&better)), BIG)) %then %do;
_hyp0="H0: pe <= pc-%trim(%left(&delta))"; _hyp1="Ha: pe > pc-%trim(%left(&delta))";
_tvalue=(_p1-_p2+&delta)/_diffSE;
_lcl = (_p1-_p2) - probit(1 - &alpha) * _diffSE;
_pvalue=(1-probnorm(_tvalue));
keep _lcl;
%end;
%else %if (%index(%quote(%upcase(&better)), SMALL)) %then %do;
_hyp0="H0: pe >= pc+%trim(%left(&delta))"; _hyp1="Ha: pe < pc+%trim(%left(&delta))";
_tvalue=(_p1-_p2-&delta)/_diffSE;
_ucl = (_p1-_p2) + probit(1 - &alpha) * _diffSE;
_pvalue=probnorm(_tvalue);
keep _ucl;
%end;
%end;
%end;
%else %if (%index(%quote(%upcase(&test)),FM)) or
(%index(%quote(%upcase(&test)),%quote(FARRINGTON)) and
%index(%quote(%upcase(&test)),%quote(MANNING))) %then %do;
%if (%quote(&delta) eq) or (%quote(&delta) eq %quote(0)) %then %do;
%let delta=0;
%end;
format _pi 11.8;
_pi = arcos(-1);
/* Note: the following is for "The smaller, the better"*/
/* Note: here _n1 and _n2 are all swapped from the original paper */
_theta1=_n2/_n1;
_a1 = 1 + _theta1;
_b1 = - (1 + _theta1 + _p1 + _theta1 * _p2 + &delta * (_theta1 + 2));
_c1 = &delta**2 + &delta * (2 * _p1 + _theta1 + 1) + _p1 + _theta1 * _p2;
_d1 = -_p1 * &delta * (1 + &delta);
_v1 = (_b1/(3*_a1))**3 - (_b1*_c1)/(6*(_a1**2)) + _d1/(2*_a1);
_u1 = sign(_v1)*sqrt((_b1/(3 * _a1))**2 - _c1/(3 * _a1));
_w1 = (1/3) * (_pi + arcos(_v1/_u1**3));
_p1_est1 = 2 * _u1 * cos(_w1) - _b1/(3 * _a1);
_p2_est1 = _p1_est1 - δ
_diffSE1=sqrt((_p1_est1 * (1 - _p1_est1) + (_p2_est1 * (1 - _p2_est1))/_theta1)/_n1);
/* Note: here _n1 and _n2 are all swapped from the original paper
for "The bigger, the better";
*/
_theta2=_n1/_n2;
_a2 = 1 + _theta2;
_b2 = - (1 + _theta2 + _p2 + _theta2 * _p1 + &delta * (_theta2 + 2));
_c2 = &delta**2 + &delta * (2 * _p2 + _theta2 + 1) + _p2 + _theta2 * _p1;
_d2 = -_p2 * &delta * (1 + &delta);
_v2 = (_b2/(3*_a2))**3 - (_b2*_c2)/(6*(_a2**2)) + _d2/(2*_a2);
_u2 = sign(_v2)*sqrt((_b2/(3 * _a2))**2 - _c2/(3 * _a2));
_w2 = (1/3) * (_pi + arcos(_v2/_u2**3));
_p2_est2 = 2 * _u2 * cos(_w2) - _b2/(3 * _a2);
_p1_est2 = _p2_est2 - δ
_diffSE2=sqrt((_p2_est2 * (1 - _p2_est2) + (_p1_est2 * (1 - _p1_est2))/_theta2)/_n2);
%if (%index(%quote(%upcase(&tail)),TWO)) or (%quote(&tail) eq 2) %then %do;
_hyp0="H0: |pc-pe|<=%trim(%left(&delta))"; _hyp1="Ha: |pc-pe| > %trim(%left(&delta))";
if (1-probnorm(abs((_p1-_p2+&delta)/_diffSE2)))>=(1-probnorm(abs((_p1-_p2-&delta)/_diffSE2))) then do;
_diffSE=_diffSE2;
_tvalue = (_p1-_p2+&delta)/_diffSE;
end;
else do;
_diffSE=_diffSE1; put _diffse1=;
_tvalue=(_p1-_p2-&delta)/_diffSE;
end;
_pvalue=(1-probnorm(abs((_p1-_p2+&delta)/_diffSE2)))+(1-probnorm(abs((_p1-_p2-&delta)/_diffSE1)));
_lcl = (_p1-_p2) - probit(1 - &alpha/2) * _diffSE;
_ucl = (_p1-_p2) + probit(1 - &alpha/2) * _diffSE;
keep _lcl _ucl;
%end;
%else %if (%index(%quote(%upcase(&tail)),ONE)) or (%quote(&tail) eq 1) %then %do;
%if (%index(%quote(%upcase(&better)), BIG)) %then %do;
_hyp0="H0: pe <= pc-%trim(%left(&delta))"; _hyp1="Ha: pe > pc-%trim(%left(&delta))";
_diffSE=_diffSE2;
_tvalue = (_p1-_p2+&delta)/_diffSE;
_lcl = (_p1-_p2) - probit(1 - &alpha) * _diffSE;
_pvalue=(1-probnorm(_tvalue));
keep _lcl;
%end;
%else %if (%index(%quote(%upcase(&better)), SMALL)) %then %do;
_hyp0="H0: pe >= pc+%trim(%left(&delta))"; _hyp1="Ha: pe < pc+%trim(%left(&delta))";
_diffSE=_diffSE1;
_tvalue = (_p1-_p2-&delta)/_diffSE1;
_ucl = (_p1-_p2) + probit(1 - &alpha) * _diffSE;
_pvalue=probnorm(_tvalue);
keep _ucl;
%end;
%end;
%end;
%else %do;
%put ==> Alert! I don%str(%')t recognize "&test"! I can only run "HA" (stand for "HAUCK ANDERSON") or;
%put +++ "FM" (stand for "FARRINGTON MANNING"). Please choose one test!;
%goto finish;
%end;
keep _p _tvalue _diffSE _hyp0 _hyp1 _pvalue;
label _tvalue="Test Statistics" _diffSE="Standard Error of the Proportion Difference"
_hyp0="Null Hypothesis" _hyp1="Alternative Hypothesis" _p="Overall Success Rate"
_n1="Sample size of Control Group" _p1="Success Rate in Control Group"
_n2="Sample size of Experimental Group" _p2="Success Rate in Experimental Group"
_s1="Number of Successes in Control Group" _s2="Number of Successes in Experimental Group"
_pvalue="P-Value" _lcl="Lower CI for (pc-pe)" _ucl="Upper CI for (pc-pe)";
run;
data %if (%quote(&outdata) ne) %then %do; &outdata %end; %else %do; %let syslast=&_tmplast_; _null_ %end;;
file print;
set _tmp1_;
format default=11.8;
%if (%quote(&test) eq) or (%quote(%upcase(&test)) eq %quote(BW)) or
(%index(%quote(%upcase(&test)),%quote(BLACKWELDER))) or
(%index(%quote(%upcase(&test)),HA)) or
(%index(%quote(%upcase(&test)),%quote(HAUCK)) and
%index(%quote(%upcase(&test)),%quote(ANDERSON))) %then %do;
put /;
%if (%index(%quote(%upcase(&tail)),TWO)) or (%quote(&tail) eq 2) %then %do;
%if (%index(%quote(%upcase(&test)),HA)) or (%index(%quote(%upcase(&test)),%quote(HAUCK)) and
%index(%quote(%upcase(&test)),%quote(ANDERSON))) %then %do;
put /@5 "Two-Sided Hypothesis Test on Two Proportions Using Normal Approximation with";
put @5 "Hauck-Anderson Correction";
put @5 76*"="/;
put @5 "Reference: Hauck w. w. and Andersen S.A. (1986) A comparison of large-sample"/
@5 "confidence interval methods for the difference of two binomial probabilities."/
@5 "AM Stat. VOL 40, 318-322.";
%end;
%else %do;
put /@5 "Two-Sided Hypothesis Test on Two Proportions Using Normal Approximation with";
put @5 "Blackwelder Procedure";
put @5 76*"="/;
put @5 "Reference: Blackwelder, W. C. (1982). Proving the null hypothesis in clinical"/
@5 "trials. Controlled Clinical Trials 3, 345-353.;";
%end;
%end;
%else %if ((%index(%quote(%upcase(&tail)),ONE)) or (%quote(&tail) eq 1)) %then %do;
put /@5 "One-Sided Hypothesis Test on Two Proportions Using Normal Approximation with";
put @5 "Blackwelder Procedure";
put @5 76*"="/;
put @5 "Reference: Blackwelder, W. C. (1982). Proving the null hypothesis in clinical"/
@5 "trials. Controlled Clinical Trials 3, 345-353.;";
%end;
%end;
%else %if (%index(%quote(%upcase(&test)),FM)) or (%index(%quote(%upcase(&test)),%quote(FARRINGTON MANNING)))%then %do;
%if (%index(%quote(%upcase(&tail)),TWO)) or (%quote(&tail) eq 2) %then %do;
put /@5 "Two-Sided Hypothesis Test on Two Proportions using Farrington Manning Procedure";
put @5 79*"="/;
%end;
%else %if ((%index(%quote(%upcase(&tail)),ONE)) or (%quote(&tail) eq 1)) %then %do;
put /@5 "One-Sided Hypothesis Test on Two Proportions using Farrington Manning Procedure";
put @5 79*"="/;
%end;
put @5 "Reference: Farrington, C. P. and Manning, G. (1990) Test statistics and sample"/
@5 "size formulae for comparative binomial trials with null typothesis of non-zero"/
@5 "risk difference or non-unity relative risk. Statistics in Medicine, VOL 9," /
@5 "1447-1454.";
%end;
put /;
put &_linesize_*"-";
%if (%index(%quote(%upcase(&tail)),TWO)) or (%quote(&tail) eq 2) %then %do;
put @5 "Hypothesis Test (Two Sided Test): ";
put @10 _hyp0 "vs " _hyp1 +(-1) ";"/;
%end;
%else %if ((%index(%quote(%upcase(&tail)),ONE)) or (%quote(&tail) eq 1)) and (%index(%quote(%upcase(&better)), BIG)) %then %do;
put @5 "Hypothesis Test (The Bigger, The Better!): ";
put @10 _hyp0 "vs " _hyp1 +(-1) ";"/;
%end;
%else %if ((%index(%quote(%upcase(&tail)),ONE)) or (%quote(&tail) eq 1)) and (%index(%quote(%upcase(&better)), SMALL))%then %do;
put @5 "Hypothesis Test (The Smaller, The Better!): ";
put @10 _hyp0 "vs " _hyp1 +(-1) ";"/;
%end;
put @5 "Overall Success Rate: " @35 _p 20.6 "; ";
put @5 "Standard Error (Diff): " @35 _diffSE 20.6 "; ";
put @5 "Test Statistics: " @35 _tvalue 20.6 "; ";
put @5 "P Value: " @35 _pvalue 20.6 "; "/;
put &_linesize_*"-"/;
%if (%index(%quote(%upcase(&tail)),TWO)) or (%quote(&tail) eq 2) %then %do;
put @5 "Two Sided Confidence Intervals:"/;
put @5 "%sysevalf((1-&alpha)*100)% CI for (pe-pc): (" _lcl +(-1) ", " _ucl +(-1) ").";
%end;
%else %if ((%index(%quote(%upcase(&tail)),ONE)) or (%quote(&tail) eq 1)) %then %do;
put @5 "One Sided Confidence Interval (Please choose one):"/;
%if (%index(%quote(%upcase(&better)), BIG)) %then %do;
put @5 "%sysevalf((1-&alpha)*100)% CI for (pe-pc): (" _lcl +(-1) ", +Inf).";
%end;
%else %if (%index(%quote(%upcase(&better)), SMALL)) %then %do;
put @5 "%sysevalf((1-&alpha)*100)% CI for (pe-pc): (-Inf, " _ucl +(-1) ").";
%end;
%end;
run;
proc datasets library=work nolist;
delete _tmp1_;
run;quit;
%finish:
run; quit;
%mend bnmleqv;
